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A method to treat a iV-component percolation model as 
effective one component model is presented by introducing a 
scaled control variable p+ . In Monte Carlo simulations on 16 3 , 
32 3 , 64 3 and 128 3 simple cubic lattices the percolation thresh- 
old in terms of p+ is determined for N = 2. Phase transitions 
are reported in two limits for the bond existence probabili- 
ties p= and p^. In the same limits, empirical formulas for 
the percolation threshold p c + as function of one component- 
concentration, fb, are proposed. In the limit p= = a new 
site percolation threshold, ffj ~ 0.145, is reported. 

PACS numbers: 64.60.-i, 64.60. Ak, 02.70.Lq, 05.10.Ln 

The percolation model goes back to Flory [EJ who in- 
troduced it in the context of polymer gelation. Since 
then it has been used in a wide range of approaches and 
techniques In standard percolation models either 

bond or site percolation is dealt with |l0|,[ll]]. Site-bond 
percolation jl2],[l3) combines the two formulations, deal- 
ing with randomly occupied sites (vertices) and randomly 
existing bonds (open edges) connecting these sites. How- 
ever in this version of the model only one active compo- 
nent exists, the other sites are considered unoccupied. A 
further generalization is to consider several components, 
which was done for site percolation as well as bond per- 
colation by Zallen [13] and called polychromatic perco- 
lation. Zallen focused on the coexistence of percolating 
species in highly connected lattices, giving a criterion for 
the occurrence of a panchromatic regime where all species 
percolate. Site-bond percolation using two components 
was investigated previously by one of us jL5| and applied 
to the question of the nuclear liquid gas phase transi- 
tion. Site-bond percolation with several species was con- 
sidered in []l6f and an approximate percolation criterion 
was given. 

In this Letter, we investigate a two component site- 
bond percolation model on a simple cubic lattice, focus- 
ing on two specific limits which exhibit novel behavior. 
Let us begin by describing the approach we have taken, 
in the general case of N different component flavors. No 
assumption concerning topological dimensions or lattice 
structure is made. We have N component concentrations 
fi with Yli=i fi = 1 an d different bond probabilities to 
connect all possible combinations of sites, resulting in 
Apar = (N — 1) + ( N+ 2 ~ ) free parameters a*. The bonds 
have been assumed to be directionless, meaning that their 
probabilities only depend on the species of the sites they 
are connecting. We now want to know in which region 



of this A pQr -dimensional parameter space an infinite net- 
work Coo of connected bonds occurs, that is, where the 
probability for a given site to belong to the infinite net- 
work, poo({ai}), is non-zero. The particular type of a 
bond shall be irrelevant in order for it to belong to the 
infinite network. For a system with N > 3 components, 
however, this approach is quite impractical. It would be 
preferable to be able to reduce the dependence of the 
order parameter, p^, to one variable at fixed particle 
concentrations. We propose one such variable in follow- 
ing the definitions of ]l5| , generalizing them to N com- 
ponents. In analogy to the bond existence probability 
used in ordinary one component bond percolation mod- 
els, which gives the bond density in the system under 
observation, we introduce the scaled control parameter 
P+- 

N 

P+= a ijPiJ- C 1 ) 

»>J=1 

Here, the pij 's denote the probability for a bond to exist 
between two sites occupied by species i and j and ay = 
otji is the probability that any given nearest neighbor 
edge is one that connects two sites of flavors i and j: 

otij = 2/,/j, (2) 

with the constraint that 53i>j a y = 1 ■ 

For our simulation we consider a two component sys- 
tem on 16 3 , 32 3 , 64 3 and 128 3 simple cubic lattices. For 
simplicity we shall call one species blue, the other red. 
We now have four free parameters to vary: The frac- 
tion of one of the components, say of the blue sites, fb 
and three bond activation probabilities: pbb for bonds 
connecting two blue sites, p rr for bonds connecting two 
red sites and p^ for 6-r-bonds. However we shall set 
p = = pbb = Prr, introducing a symmetry in the sys- 
tem. This is motivated by considerations of, for example, 
isospin symmetry, where the e + e + and e~e~ interactions 
are identical. Eqs. ([!]) and (|2|) now read 

P+ =a^p^ + a= p=, (3) 

with a = = (1 — a^), and 

= 2f b (l - f b ), (4) 

respectively, where we have replaced the double indices 
by a more intuitive notation for only two components. 
Again the question is in which region of the three di- 
mensional p— — p^ — fb space an infinite network of bonds 
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appears in the lattice. We defer the question of the con- 
centration dependence to later and for the moment set fb 
to some fixed value, which shall, for now, be fb = 0.5. In 
the simulation the lattice is populated at random, with- 
out correlations, according to fb and bonds are formed for 
varied values of p =1 p^ £ [0, 1] using a Monte Carlo algo- 
rithm. The resulting cluster structure is analyzed using 
a cluster- find-algorithm described in || and Poc(p=,P^) 
is recorded. As always, a cluster is defined as a set of 
vertices connected by open edges. In Fig. [I] we show 
as a function of the two control parameters, p = and p^, 
at a concentration of fb — 0.5. We can see that p^ 
changes from (front corner) to 1 (back corner), with 
a critical line of a second order phase transition in the 
{p=,P^) -plane. We now follow our previous considera- 
tion and analyze the same data in terms of the scaled 
control parameter p + , the result of which is displayed 
in Fig. g. Two distinct branches are seen, both in the 
shape of a second order phase transition. Figure || first 
suggests that p+ is a good control parameter, as we only 
have these two universal scaled curves, but furthermore 
that another transition seems to take place in the system. 
Analysis of the data shows that the 'upper' branch con- 
stitutes of points with both p = and p^ non-zero, whereas 
all points with (p = = 0,p^ ^ 0) and (p = ^ 0,p=i = 0) 
fall on the 'lower' branch. For other values of fb the same 
behavior is found, but then the two curves for the zero- 
limits in p = or p^ are not the same. Also shown in Fig. 
^| are the expectations from one component bond per- 
colation theory, p^ oc {p + — p c + Y , (smooth lines). For 
the finite (p = ,p^) regime the critical value p c + has the 
same numerical value as the bond existence probability 
in one component bond percolation (aside from a small 
difference due to finite size effects), p\ = 0.251 ± 0.002. 
As long as both bond-types are active, the system under 
observation here and the one component model show an 
identical phase transition behavior, which is consistent 
with the findings presented in Gaj. In the zero-limits of 
p= or p^, however, p c + is shifted to p c + = 0.280 ± 0.002. 
The critical exponent /3 = 0.41 from one component bond 
percolation theory is the same in both cases presented 
here, as shown in the double logarithmic plot in the inset 
of Fig. |. 

What causes the change of p c + in the limits p = — > 
and p^ — > 0? We restrict ourselves to a discussion of 
p^ — > 0, as the same line of arguments applies in the 
other limit. For simplicity we first choose p— = 1. Then 
setting p^L = corresponds to a lattice in which all avail- 
able bb- and rr-edges are open, but all &r-edges are closed. 
There will be an infinite network present, as usual iden- 
tified with the biggest percolating cluster (with two com- 
ponents each above the site percolation threshold being 
present we could have two percolating clusters), with p^ 
given by 



^ = o) = f b -J2 



< s, 



(5) 



s=l 



where n b s is the number of ^-clusters of size s and the 
upper limit in the sum is meant to indicate that the in- 
finite cluster is excluded. Letting p^ — e = 1/A e£ j ses 
introduces on average one 6r-bond in the lattice. This 
bond might, with some small probability, which is re- 
lated to the number of perimeter edges of all finite clus- 
ters, connect two finite clusters. With a noticeably higher 
probability however, it will open an edge that connects 
the infinite cluster to the biggest r-cluster, C r max , with 
\C r max \ ~ | Coo | . The infinite cluster can now consist of 
both b- and r-sites and Poa reads 



Poc(p# 



oc — 1 

1 ~~ ^2 Us s ' 



(6) 



which, by rewriting /& in Eq. |^ as fb = 1 — J2T=i n s s ' 
leads to a difference in j>oo : 

doc (e) = Poo (fV = e ) - Poo (p=£ = 0) 



J2 < s + J2 < s) - ]T n s s, (7) 



with the property that lim c _»o (e) ^ for all events 
where Coo and are connected. Averaging over all 

events will yield some effective 5^ with < < Soo- 
This leads to the conclusion that in the limit p^ — » a 
first order phase transition takes place in j»oo ■ The sim- 
ulation results support this conjecture. Figure [| shows 
these results, where the average was taken over events 
in which the number of 6r-bonds actually formed in the 
simulation, Ufe r , was non-zero. For values of p— other 
than 1 one still finds the same behavior, somewhat less 
pronounced due to smaller |Coo| and \C'! max \. The same 
holds for fb 7^ 0.5 where it is clear that the effect vanishes 
continuously in the limit of a one component bond perco- 
lation system: fb — > or — > 1. Obviously the question 
arises what the behavior of p\ is. We found that for 
p^ — > the percolation threshold p c + continually changes 
to a new, concentration dependent value. In Fig. ||(a) 
we show this transition, which is of second order type. 
We also show a fit to the curve, which is of the form 



p c + {p=,p^ -> 0) 



(u + v p-l) 



+ t, 



(8) 



where the fit parameters were found to be u — 34.6 ±0.3, 
v = 1823 ± 72 and t = 0.246 ± 0.001. It has to be noted 
that this formula only stands on empirical grounds, fits 
by exponential functions may also be useful. As stated 
earlier, the same line of thought is applicable to the p = — > 
limit and indeed do our simulations present the same 
results in the system with fb = 0.5. For fb ^ 0.5 this 
direct symmetry is broken, but qualitatively the results 
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are still the same. Even for fb = 0.5 they differ however 
with respect to the afore-mentioned dependence of p c + on 
the concentration fb in the limits p = = and = 0. 
For p+(fb,Pj£ — 0) we find the functional form depicted 
in Fig. ^(b). The results shown are fits of p L + to the 
scaling relation, 



site percolation model which yields /f; 



0.3116. In 



\p\(L)-p%\^L- 



-1/u 



Fig. ||(b) determines the critical point of the phase 
transition line p+(fb,P= — 0). Another argument in fa- 
vor of a qualitatively new behavior arises by considering 
the density of accessible edges at the critical concentra- 
tions in the two models. For the one considered here it is 
given by a 7 t(f b ), see Eq. (|J), for the usual site percola- 
(9) tion model we might define one in an analogous manner: 



as given in Ref. po| |, where we kept v fixed at 0.88 and 
lattice sizes L — 16, L = 32, L = 64 and L = 128 were 
taken into account. The simulation data is fitted with 



p c +(fb,P? = 0) 



2(/ 6 



h + mfb 



(10) 



This is Eq. (|J) with the purely empirical assumption 
of a hyperbola for p'L(fb), which is in agreement with 
the results of Heermann and Stauffer for a one compo- 
nent site-bond model 13 1. Fitting the parameters to our 
simulation data results in h = 4.007 ± 0.002 and m — 
—4.428 ± 0.005. A comparison with the formula given 
by Heermann and Stauffer yields h = l/.Pb ond (0) = 4.019 
and m = (1 - P c bond ) / (jp%(f c site - 1)) = -4.386, where 
Pbond = 0.2488 and f^ ite = 0.3116 are the percolation 
thresholds for one component bond- and site- percolation 
on a three dimensional simple cubic lattice, respectively. 
However, in contrast to the one component model, due 
to the symmetry introduced in the system, notably that 
p = = pbb = Prr, we only have f b G [0.0, 0.5] as indepen- 
dent regime here, with the interval fb G [0.5, 1.0] being 
symmetric to the one shown here, with only the roles of 
blue and red sites being switched. This is a manifesta- 
tion of the two components behaving like two superposed, 
non-interfering one component site-bond percolation sys- 
tems and lets us conclude that the phase transition in p c + 
for p^ — > discussed above may be interpreted as an ef- 
fective transition from a one component bond percolation 
model to a one component site-bond percolation model. 
For the case of the second limit discussed, p = — > 0, the 
results are shown in Fig. ||. In Fig. ||(a) we go back to 
the bond existence probability p^ . Its critical value, as a 
function of the concentration, p^ifb), is well reproduced 
with an exponential fit: 



V = P^t(fb) = a exp(-d fb) + • 



(11) 



Fitting this empirical formula to the simulation data 
gives a = 2.1 ±0.07, d = 10.9 ±0.2 and c = 0.547 ±0.002. 
We define f? = (p^ 1 (p^^l) and numerically get f£ ~ 0.14 



from Eq. (11), whereas in an independent simulation we 
find fi = 0.145 ± 0.001, which again is a result from a 



fit to the scaling relation Eq. 



The parameter f? can 



be regarded as a new threshold, which would correspond 
to the percolation threshold in a simple site percolation 
model in which nearest neighbors only belong to the same 
cluster if they are of opposite flavor, unlike the normal 



(f c site y. By setting a # (/ 6 °) 



O-sitt 



one would 



expect f§ to be 0.051 which stands in contradiction to 
our findings. 

In concluding, we introduced a new way to treat N- 
component percolation. This approach was applied to a 
two component site-bond percolation model and new first 
order phase transitions of p^ were reported in the limits 
p^ — > and p = — » 0. In the latter case we could further- 
more establish a novel empirical formula for the percola- 
tion threshold as a function of component concentration, 
whereas in the first earlier findings of a one component 
site-bond percolation model were found to apply in the 
two component model too. The field for future work 
in this area seems vast, one might, for example, try to 
apply the same method to multi-component systems on 
lattices of higher dimensions and/or higher connectivity. 
This approach should also find a broad range of possi- 
ble applications. One might think of special networks or 
gelation phenomena with several components involved, 
which only interact with each other, as well as wetting 
phenomena. Furthermore an application to stock-market 
simulations seems possible and is being undertaken by 
the authors. 
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FIG. 1. Probability to belong to the infinite cluster, p^, 
in a simple cubic two component site-bond percolation lattice 
of size 128 3 as a function of the two parameters, p= and p^, 
calculated for a fraction of the blue species of ft = 0.5. 




0.8 



FIG. 2. Probability to belong to the infinite cluster, p^, 
as a function of the scaled control parameter, p+, in a 128 3 
simple cubic lattice with ft = 0.5, for values of p=,p_^ £ [0, 1], 
fitted with poo oc (p+ — P+) 13 (solid line). Here p\ ~ 0.251 
for the 'upper' branch, p c + ~ 0.280 for the 'lower' branch and 
(3 = 0.41 in both cases. The inset shows the same data in a 
double logarithmic representation. 




1 2 3 _ 4 5 

A = (p+- 0.5) 10 7 
FIG. 3. Probability to belong to the infinite cluster, p^ , as 
a function of the scaled control parameter p + , in a 128 3 simple 
cubic lattice with fb = 0.5, at fixed p= = 1 for varied p^. 50 
independent simulation events where taken into account and 
the average was taken over events in which the number of 
6r-bonds actually formed in the simulation, nj, r , was non-zero. 
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P* fb 
FIG. 4. Critical value of the scaled control parameter, p+, 

in a 128 3 simple cubic lattice with ft = 0.5, (a) plotted as a 

function of p^, and (b) plotted as a function of the fraction 

of blue sites, /j,, in the limit p^ = and obtained by a fit to 

the scaling law |p+(L) — p+| oc L^ 1 ^. The errors, estimated 

as described in H, are smaller than the symbol sizes in (b). 
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FIG. 5. Critical value of p^ and of the scaled control pa- 
rameter, and p c + , plotted as a function of the fraction of 
blue sites /;, in the limit p^ — 0. The data has been obtained 
by a fit to \p+(L) — p c + \ oc L~ 1//,y . Again the errors, estimated 
as described in M, are smaller than the symbol sizes. 



5 



